
#################################################################################################

#  Description: Generating Tables and Figures from Study 2 in Fleming, Thomas G. (Forthcoming) 
#               'Partisanship and the Effectiveness of Personal Vote-Seeking', Legislative Studies 
#               Quarterly.
#  Author:      Thomas G. Fleming
#  Date:        18/01/2021
#  Note:        Requires files "study2_data.csv" to be saved in same folder

#################################################################################################

rm(list=ls())

library(stargazer)
library(coefplot)
library(gridExtra)  # for grid.arrange()
library(grid)       # for textGrob()
library(arsenal)    # for tableby()
library(xtable)

library(rstudioapi)
current_path <- getActiveDocumentContext()$path
setwd(dirname(current_path))

######################################################
# Loading data and setting reference category of IV  #
######################################################

  data <- read.csv("./study2_data.csv")
  data$treatment <- relevel(data$treatment, ref = "Control")

##############################################################
# Figure 4 - Voter responses to MP's activity (full sample)  #
##############################################################

  # running main model
  
  ols.fullsample <- lm(vote ~
                         treatment,
                       data = data) 
  
  # creating coefficient plot from this model
  
  coefplot(ols.fullsample,
           title = NULL,
           xlab = "\n Average treatment effects",
           ylab = NULL,
           innerCI = 1.96, # 95% CI
           lwdInner = 0.5,
           outerCI = 0,
           intercept = FALSE,
           pointSize = 2,
           color = "black",
           zeroColor="black",
           zeroLWD=0.5,
           newNames = c("treatmentInactive MP" = "Inactive MP\n treatment",
                        "treatmentActive MP" = "Active MP\n treatment"))+
    theme_classic()+
    theme(axis.title.x = element_text(face = "bold"),
          panel.border = element_rect(colour = "black", fill = NA, size = 0.7), # adds rectangle
          axis.line = element_blank())                                          # removes axis lines because they reinforce the rectangle 
  ggsave("./figure4.jpg", width = 6, height = 2.5, units = "in")
  
##############################################################
# Table C4 - OLS regression of support for MP (full sample)  #
##############################################################

  # creating regression table from above model
  
  stargazer(ols.fullsample,
            type = "html",
            title = "Table C4: OLS regression of support for MP (full sample)",
            dep.var.caption = "Support for MP",
            dep.var.labels.include = FALSE,
            model.numbers = FALSE,
            column.labels = c("(C3)"),
            out = "./tableC4.html",
            covariate.labels = c("Active MP treatment", "Inactive MP treatment", "Constant"),
            omit.stat = c("f", "ser", "adj.rsq"),
            no.space = TRUE)
  
#################################################################
# Figure 5 - Voter responses to MP's activity (by partisanship) #
#################################################################

  # running model separately by two partisanship groups
  
  nonpartisan.ols <- lm(vote ~
                          treatment,
                        data = data[(data$partisanship.strength == "Not very strong" | 
                                     data$partisanship.strength == "None"),])                                
  
  partisan.ols <- lm(vote ~
                       treatment,
                     data = data[(data$partisanship.strength == "Very strong" |
                                  data$partisanship.strength == "Fairly strong"),])    
  
  # creating and combining coefficient plots from these models
  
  inactive.twogroups <- multiplot(partisan.ols,
                                  nonpartisan.ols,
                                  coefficients = c("treatmentInactive MP"),
                                  title = "(a) Inactive MP treatment",
                                  xlab = NULL,
                                  ylab = NULL,
                                  innerCI = 1.96, # 95% CI
                                  lwdInner = 0.5,
                                  outerCI = 0,
                                  intercept = FALSE,
                                  pointSize = 2,
                                  zeroColor="black",
                                  zeroLWD=0.5,
                                  names = c("Partisan\n(n = 663)", "Non-partisan\n(n = 685)"),
                                  by = "Model")+
    theme_classic()+
    theme(plot.title = element_text(size = 10),
          legend.position = "none",
          panel.border = element_rect(colour = "black", fill = NA, size = 0.7), # adds rectangle
          axis.line = element_blank())+                                         # removes axis lines because they reinforce the rectangle
    scale_color_manual(values = c("black", "black"))
  
  active.twogroups <- multiplot(partisan.ols,
                                nonpartisan.ols,
                                coefficients = c("treatmentActive MP"),
                                title = "(b) Active MP treatment",
                                xlab = NULL,
                                ylab = NULL,
                                innerCI = 1.96, # 95% CI
                                lwdInner = 0.5,
                                outerCI = 0,
                                intercept = FALSE,
                                pointSize = 2,
                                zeroColor="black",
                                zeroLWD=0.5,
                                names = c("Partisan\n(n = 663)", "Non-partisan\n(n = 685)"),
                                by = "Model")+
    theme_classic()+
    theme(plot.title = element_text(size = 10),
          legend.position = "none",
          panel.border = element_rect(colour = "black", fill = NA, size = 0.7), # adds rectangle
          axis.line = element_blank())+                                         # removes axis lines because they reinforce the rectangle
    scale_color_manual(values = c("black", "black"))
  
  figure5 <- grid.arrange(inactive.twogroups, 
                          active.twogroups, 
                          ncol = 2, 
                          bottom = textGrob("\n               Conditional average treatment effects",
                                            gp = gpar(fontface = "bold",
                                                      fontsize = 11))) 
  ggsave("./figure5.jpg", figure5, width = 6, height = 3.5, units = "in")
  
###############################################################
# Table C5 - OLS regression of support for MP (two-way split) #
###############################################################
  
  # repeating above models with separate (but identical) DVs for easier table formatting
  
  data$vote1 <- data$vote
  
  nonpartisan.ols.fortable <- lm(vote1 ~
                                   treatment,
                                 data = data[(data$partisanship.strength == "Not very strong" | 
                                              data$partisanship.strength == "None"),])                                
  
  partisan.ols.fortable <- lm(vote ~
                                treatment,
                              data = data[(data$partisanship.strength == "Very strong" |
                                           data$partisanship.strength == "Fairly strong"),])  
  
  data$vote1 <- NULL
  
  # creating regression table from above models
  
  stargazer(partisan.ols.fortable,
            nonpartisan.ols.fortable,
            type = "html",
            title = "Table C5: OLS regression of support for MP (two-way split)",
            dep.var.caption = "Support for MP",
            dep.var.labels.include = TRUE,
            dep.var.labels = c("Partisan", "Non-partisan"),
            model.numbers = FALSE,
            column.labels = c("(C4)", "(C5)"),
            out = "./tableC5.html",
            covariate.labels = c("Active MP treatment", "Inactive MP treatment", "Constant"),
            omit.stat = c("f", "ser", "adj.rsq"),
            no.space = TRUE)
  
#################################################################
# Figure 6 - Voter responses to MP's activity (by partisanship) #
#################################################################

  # running model separately by four partisanship groups
  
  D.noparty.ols <- lm(vote ~
                        treatment,
                      data = data[(data$partisanship.strength == "None"),]) 
  
  C.notverystrong.ols <- lm(vote ~
                              treatment,
                            data = data[(data$partisanship.strength == "Not very strong"),])  
  
  B.fairlystrong.ols <- lm(vote ~
                             treatment,
                           data = data[(data$partisanship.strength == "Fairly strong"),]) 
  
  A.verystrong.ols <- lm(vote ~
                           treatment,
                         data = data[(data$partisanship.strength == "Very strong"),])    
  
  # creating and combining coefficient plots from these models
  
  inactive.fourgroups <- multiplot(D.noparty.ols, 
                                   C.notverystrong.ols, 
                                   B.fairlystrong.ols, 
                                   A.verystrong.ols,
                                   coefficients = c("treatmentInactive MP"),
                                   title = "(a) Inactive MP treatment",
                                   xlab = NULL,
                                   ylab = NULL,
                                   innerCI = 1.96,  # 95% CI
                                   lwdInner = 0.5,
                                   outerCI = 0,
                                   intercept = FALSE,
                                   pointSize = 2,
                                   zeroColor = "black",
                                   zeroLWD = 0.5,
                                   names = c("  None\n(n = 310)",
                                             " Not very strong\n(n = 375)",
                                             "Fairly strong\n(n = 505)",
                                             "Very strong\n(n = 158)"),
                                   by = "Model")+
    theme_classic()+
    theme(plot.title = element_text(size = 10),
          #  axis.text.y = element_text(color = "white"),
          legend.position = "none",
          panel.border = element_rect(colour = "black", fill = NA, size = 0.7), # adds rectangle
          axis.line = element_blank())+                                         # removes axis lines because they reinforce the rectangle
    scale_color_manual(values = c("black", "black", "black", "black"))
  
  active.fourgroups <- multiplot(D.noparty.ols, 
                                 C.notverystrong.ols, 
                                 B.fairlystrong.ols, 
                                 A.verystrong.ols,
                                 coefficients = c("treatmentActive MP"),
                                 title = "(b) Active MP treatment",
                                 xlab = NULL,
                                 ylab = NULL,
                                 innerCI = 1.96,  # 95% CI
                                 lwdInner = 0.5,
                                 outerCI = 0,
                                 intercept = FALSE,
                                 pointSize = 2,
                                 zeroColor = "black",
                                 zeroLWD = 0.5,
                                 names = c("  None\n(n = 310)",
                                           " Not very strong\n(n = 375)",
                                           "Fairly strong\n(n = 505)",
                                           "Very strong\n(n = 158)"),
                                 by = "Model")+
    theme_classic()+
    theme(plot.title = element_text(size = 10),
          legend.position = "none",
          panel.border = element_rect(colour = "black", fill = NA, size = 0.7), # adds rectangle
          axis.line = element_blank())+                                         # removes axis lines because they reinforce the rectangle
    scale_color_manual(values = c("black", "black", "black", "black"))
  
  figure6 <- grid.arrange(inactive.fourgroups, 
                          active.fourgroups, 
                          ncol = 2, 
                          bottom = textGrob("\n               Conditional average treatment effects",
                                            gp = gpar(fontface = "bold",
                                                      fontsize = 11))) 
  ggsave("./figure6.jpg", figure6, width = 6, height = 4, units = "in")
  
################################################################
# Table C6 - OLS regression of support for MP (four-way split) #
################################################################
  
  # repeating above models with separate (but identical) DVs for easier table formatting
  
  data$vote1 <- data$vote
  data$vote2 <- data$vote
  data$vote3 <- data$vote
  
  D.noparty.ols.fortable <- lm(vote1 ~
                                 treatment,
                               data = data[(data$partisanship.strength == "None"),])  
  
  C.notverystrong.ols.fortable <- lm(vote2 ~
                                       treatment,
                                     data = data[(data$partisanship.strength == "Not very strong"),])  
  
  B.fairlystrong.ols.fortable <- lm(vote3 ~
                                      treatment,
                                    data = data[(data$partisanship.strength == "Fairly strong"),]) 
  
  A.verystrong.ols.fortable <- lm(vote ~
                                    treatment,
                                  data = data[(data$partisanship.strength == "Very strong"),])   
 
  data$vote1 <- NULL
  data$vote2 <- NULL
  data$vote3 <- NULL
  
  # creating regression table from above models
  
  stargazer(A.verystrong.ols.fortable,
            B.fairlystrong.ols.fortable,
            C.notverystrong.ols.fortable,
            D.noparty.ols.fortable,
            type = "html",
            title = "OLS regression of support for MP (four-way split)",
            dep.var.caption = "Support for MP",
            dep.var.labels.include = TRUE,
            dep.var.labels = c("Very strong id", "Fairly strong id", 
                               "Not very strong id", "No id"),
            model.numbers = FALSE,
            column.labels = c("(C6)", "(C7)", "(C8)", "(C9)"),
            out = "./tableC6.html",
            covariate.labels = c("Active MP treatment", "Inactive MP treatment", "Constant"),
            omit.stat = c("f", "ser", "adj.rsq"),
            no.space = TRUE)
  
####################################################
# Table C1 - Summary statistics by treatment group #
####################################################
  
  # creating summary table
  
  randomisation.table <- tableby(treatment ~ 
                                   age + 
                                   gender +
                                   socialgrade +
                                   GEvote2017 +
                                   refvote2016 +
                                   partisanship +
                                   partisanship.strength +
                                   political.attention, 
                                 data = data, 
                                 numeric.stats=c("mean"),   
                                 test = FALSE,            # turns off test column
                                 total = FALSE)           # turns off total column
  randomisation <- xtable(as.data.frame(summary(randomisation.table, text = TRUE)))
  
  # formatting summary table
  
  caption(randomisation) <- "Table C1: Summary statistics by treatment group"
  
  randomisation[,1] <- c("", "Age (mean)", 
                         "Gender", "--- Female", "--- Male", 
                         "Social grade", "--- A/B", "--- C1/C2", "--- D/E",
                         "2017 GE vote ", "--- Conservative", "--- Didn't vote", "--- Labour", "--- Other (incl. DK)",
                         "2016 referendum vote", "", "--- Didn't vote", "--- Leave", "--- Remain",
                         "Party id", "--- Conservative", "--- Labour", "--- None", "--- Other",
                         "Party id (strength)", "--- Fairly strong", "--- None", "--- Not very strong", "--- Very strong",
                         "", "Political attention (mean)")
  
  randomisation <- randomisation[c(2,3,4,5,6,7,8,9,10,11,13,14,12,15,18,19,17,20,21,22,24,23,25,29,26,28,27,31),]
  
  # printing summary table
  
  print(randomisation,
        type = "html",
        file = "./tableC1.html",
        include.rownames = FALSE,
        size= "scriptsize",
        hline.after = c(-1,-1,0,28,28))
  
###############################################################################
# Table C2 - Proportion of responses to manipulation check by treatment group #
###############################################################################
  
  manipulation <- xtable(prop.table(table(Treatment = data$treatment, Response = data$perceived.MP.effort), 1))
  
  align(manipulation) <- "lcccc"
  caption(manipulation) <- "Table C2: Proportion of responses to manipulation check by treatment group"
  
  print(manipulation,
        type = "html",
        file = "./tableC2.html",
        hline.after = c(-1, -1, 0, 3, 3))
  
#######################################################################
# Table C3 - Partisanship's effect on support for MP (OLS regression) #
#######################################################################
  
  # setting reference category of partisanship variable
  
  data$partisanship <- relevel(data$partisanship, ref = "None")
  
  # modelling respondents' answers by their party identification
  
  labour.check1 <- lm(vote ~
                        partisanship,
                      data = data) 
  
  labour.check2 <- lm(vote ~
                        partisanship +
                        treatment,
                      data = data) 
  
  # producing regression table of the above models
  
  stargazer(labour.check1,
            labour.check2,
            type = "html",
            title = "Table C3: Partisanship's effect on support for MP (OLS regression)",
            dep.var.caption = "Support for MP",
            dep.var.labels.include = FALSE,
            model.numbers = FALSE,
            column.labels = c("(C1)", "(C2)"),
            out = "./tableC3.html",
            covariate.labels = c("Conservative id", 
                                 "Labour id", 
                                 "Other party id", 
                                 "Active MP treatment", 
                                 "Inactive MP treatment", 
                                 "Constant"),
            omit.stat = c("f", "ser", "adj.rsq"),
            no.space = TRUE)
  
##################################################################
# Table C7 - OLS regression of support for MP (with interaction) #
##################################################################
  
  # creating continuous partisanship variable for inclusion in model
  
  data$partisanship.strength.continuous <- ifelse(data$partisanship.strength == "Very strong", 3,
                                           ifelse(data$partisanship.strength == "Fairly strong", 2,
                                           ifelse(data$partisanship.strength == "Not very strong", 1,
                                           ifelse(data$partisanship.strength == "None", 0, NA))))
  
  # models treatment effects in one model, with an interaction
  
  combined.model <- lm(vote ~
                         treatment*partisanship.strength.continuous,
                       data = data)

  # creates a regression table of this model
  
  stargazer(combined.model,
            type = "html",
            title = "Table C7: OLS regression of support for MP (with interaction)",
            dep.var.caption = "Support for MP",
            dep.var.labels.include = FALSE,
            model.numbers = FALSE,
            column.labels = c("(C10)", "(C11)", "(C12)", "(C13)"),
            out = "./tableC7.html",
            covariate.labels = c("Active MP treatment", "Inactive MP treatment", "Partisanship",
                                 "Active MP treatment * Partisanship", "Inactive MP treatment * Partisanship", 
                                 "Constant"),
            omit.stat = c("f", "ser", "adj.rsq"),
            no.space = TRUE)
  
######################################################################################
# Figure C1 - Marginal effect of active/inactive MP treatments on support for the MP #
######################################################################################
  
  # plots marginal effect of treatments in above model at different levels of partisanship
  
  library(interplot)
  
  interplot(m = combined.model, var1 = "treatment", var2 = "partisanship.strength.continuous", ci = 0.95, 
            facet_labs = c("Active MP treatment", "Inactive MP treatment"))+
    geom_hline(yintercept = 0, linetype = "dashed")+
    theme_classic()+
    labs(x = "\nVoter partisanship", y = "Marginal effect", caption = "")+
    theme(axis.title.x = element_text(face = "bold"),
          axis.title.y = element_text(face = "bold"),
          panel.border = element_rect(colour = "black", fill = NA, size = 0.7),
          panel.background = element_blank(),
          axis.line = element_blank())+
    scale_x_continuous(breaks = c(0,1,2,3), labels = c("None", "Not very strong", "Fairly strong", "Very strong            "))
  dev.off()
  
  ggsave("./figureC1.jpg", width = 7, height = 4, units = "in")
  
##########################################################################################
  
  rm(list=ls())
  